use "/EDcomp.dta"
merge 1:1 country_id  year using "/nocorrupt.dta"
drop _merge
sort country_id year
merge 1:1 country_id year using "/health.dta"
drop _merge
xtset country_id year
set matsize 5000

///run the model
forvalues i = 1/900 {

*Print out an iteration number
display `i'

*Fit the model, using the ith draw from the UDS posterior
quietly xtreg F.e_peinfmor_ln EDcomp_v`i' nocorrupt_v`i' e_migdppc_ln e_miurbani e_migdpgro L.e_peinfmor_ln e_peinfmor_ln_rt i.year, fe vce (clu country_id)

*Extract the coefficients and variance-covariance matrix
matrix b = e(b)
matrix V = e(V)
local blength = colsof(b)
matrix rsq = e(r2)

*Preserve the dataset, take a single multivariate normal draw from the
*posterior distribution of the coefficients, and restore the dataset.
*We use the capture command to catch possible errors in drawnorm
*and drop these iterations gracefully.
preserve 
capture quietly drawnorm b1-b`blength', double n(1) means(b) cov(V) clear
if _rc==0 {
mkmat b1-b`blength', matrix(bsample)
matrix posterior = nullmat(posterior) \ bsample
matrix rsquared = nullmat(rsquared) \ rsq
}
else {
display "Error drawing sample...iteration dropped"
}
restore

*Closes the Monte Carlo Loop
}

*Get posterior ready to work with
svmat posterior
svmat rsquared

*Calculate means and standard deviations
tabstat posterior*, stat(mean sd)

*Find the bounds of the 95 percent credible interval
centile posterior*, centile(2.5, 97.5)

* Find the R-Squared
tabstat rsquared*, stat(mean sd)


